rm(list = ls())
library(foreign)
library(survival)

setwd("~/Dropbox/Wikipedia/final-data/Data Archive/") #Set to directory
wiki <- read.dta("Kalla_Aronow_2015.dta")

wiki$dead <- is.na(wiki$time_alive)
wiki$time_alive2 <- wiki$time_alive
wiki$time_alive2[wiki$dead] <- max(wiki$time_alive*2,na.rm=TRUE)
wiki$time_alive2[wiki$time_alive==0] <- 0.1
wiki$hours_alive <- wiki$time_alive2*2.77778e-7

wiki$survdata <- Surv(time=wiki$hours_alive)#,event=dead)

####### main plot
m <- matrix(c(1,2,3,4,4,4),nrow = 2,ncol = 3,byrow = TRUE)
layout(mat = m,heights = c(0.8,0.3,0))
par(mar = c(0,3,3,1))

plot_colors2 <- c("black", "red", "black","red")
plot_lty2 = c(1,1,2,2)
plot_colors <- c("black","black", "red", "red")
plot_lty = c(1,2,1,2)
plot_colorsC <- c("black","red")
plot_lwd=1.25

maxtime <- 250

fit1 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=(wiki$round==1 | wiki$round==4))
plot(fit1,xlim=c(0,maxtime),col=plot_colors,lty=plot_lty,lwd=plot_lwd, main="Both Cited and Uncited Edits\n Active Senators (Studies 1 and 4)")
mtext(text="Proportion Surviving", side=2, line=2, cex=0.7)
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit1$strata)) { #number of curves
     temp <- fit1[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colors[i],lty=plot_lty[i])
     }
par(mar = c(0,2,3,1))

fit2 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=wiki$round==2)
plot(fit2,xlim=c(0,maxtime),col=plot_colorsC,lty=plot_lty2,lwd=plot_lwd, main="Only Cited Edits\n Active Senators (Study 2)")
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit2$strata)) { #number of curves
     temp <- fit2[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colorsC[i],lty=plot_lty2[i])
     }
par(mar = c(0,2,3,1))

fit3 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=wiki$round==3)
plot(fit3,xlim=c(0,maxtime),col=plot_colorsC,lty=plot_lty2,lwd=plot_lwd, main="Only Cited Edits\n Retired/Deceased Senators (Study 3)")
for (i in 1:length(fit3$strata)) { #number of curves
     temp <- fit3[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colorsC[i],lty=plot_lty2[i])
     }
mtext(text="Hours", side=1, line=2,cex=0.7)
par(mar = c(0,1.5,3,1))

plot(1, type = "n", axes=FALSE, xlab="", ylab="")

legend(x = "top",inset = 0,
       legend = c("Positive Cite", "Negative Cite", "Positive Uncite", "Negative Uncite"), 
       col=plot_colors2, lty=plot_lty2, lwd=2, cex=1.2, horiz = TRUE)
       
##########
dev.off()
m <- matrix(c(1,2,3,4,5,6,7,7,7,7,7,7),nrow = 3,ncol = 3,byrow = TRUE)
layout(mat = m,heights = c(0.3,0.3,0.15,0))
par(mar = c(3,3,3,1))

# individual plots

fit1 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=(wiki$round==1))
plot(fit1,xlim=c(0,maxtime),col=plot_colors,lty=plot_lty,lwd=plot_lwd, main="Both Cited and Uncited Edits\n Active Senators (Study 1)")
mtext(text="Proportion Surviving", side=2, line=2, cex=0.7)
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit1$strata)) { #number of curves
     temp <- fit1[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colors[i],lty=plot_lty[i])
     }

fit2 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=wiki$round==2)
plot(fit2,xlim=c(0,maxtime),col=plot_colorsC,lty=plot_lty2,lwd=plot_lwd, main="Only Cited Edits\n Active Senators (Study 2)")
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit2$strata)) { #number of curves
     temp <- fit2[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colorsC[i],lty=plot_lty2[i])
     }

fit3 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=wiki$round==3)
plot(fit3,xlim=c(0,maxtime),col=plot_colorsC,lty=plot_lty2,lwd=plot_lwd, main="Only Cited Edits\n Retired/Deceased Senators (Study 3)")
for (i in 1:length(fit3$strata)) { #number of curves
     temp <- fit3[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colorsC[i],lty=plot_lty2[i])
     }
mtext(text="Hours", side=1, line=2,cex=0.7)

fit4 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=(wiki$round==4))
plot(fit4,xlim=c(0,maxtime),col=plot_colors,lty=plot_lty,lwd=plot_lwd, main="Both Cited and Uncited Edits\n Active Senators (Study 4)")
mtext(text="Proportion Surviving", side=2, line=2, cex=0.7)
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit1$strata)) { #number of curves
     temp <- fit4[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colors[i],lty=plot_lty[i])
     }

fit5 <- survfit(wiki$survdata~factor(wiki$treat_final),subset=(wiki$round==5 & wiki$treatment_post_order <= 50))
plot(fit5,xlim=c(0,maxtime),col=plot_colors,lty=plot_lty,lwd=plot_lwd, main="Both Cited and Uncited Edits\n Active Senators (Study 5, 1st Half)")
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit1$strata)) { #number of curves
     temp <- fit5[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colors[i],lty=plot_lty[i])
     }
     
fit5b <- survfit(wiki$survdata~factor(wiki$treat_final),subset=(wiki$round==5 & wiki$treatment_post_order > 50))
plot(fit5b,xlim=c(0,maxtime),col=plot_colorsC,lty=plot_lty2,lwd=plot_lwd, main="Only Cited Edits\n Active Senators (Study 5, 2nd Half)")
mtext(text="Hours", side=1, line=2,cex=0.7)
for (i in 1:length(fit2$strata)) { #number of curves
     temp <- fit5b[i]
     lines(c(max(temp$time), maxtime), rep(min(temp$surv),2),col=plot_colorsC[i],lty=plot_lty2[i])
     }

plot(1, type = "n", axes=FALSE, xlab="", ylab="")

legend(x = "top",inset = 0,
       legend = c("Positive Cite", "Negative Cite", "Positive Uncite", "Negative Uncite"), 
       col=plot_colors2, lty=plot_lty2, lwd=2, cex=1.2, horiz = TRUE)

